Not to appear in Nonlearned J., 45. 



Phosphorous chemistry in the shocked region L1157 Bl 

T. Aota and Y. Aikawa 
Department of Earth and Planetary Sciences, Kobe University, Kobe 657-8501, Japan 

ABSTRACT 

We study the evolution of phosphorous-bearing species in one-dimensional 
C-shock models. We find that the abundances of P-bearing species depend sen- 
sitively on the elemental abundance of P in the gas phase and on the abundance 
i-^ ■ of N atoms in the pre-shock gas. The observed abundance of PN and the non- 

detection of PO towards L1157 Bl are reproduced in C-shock models with shock 
velocity v = 20 km s" 1 and pre-shock density n(H 2 ) = 10 4 — 10 5 cm -3 , if the ele- 
mental abundance of P in the gas phase is ~ 10 -9 and the N-atom abundance is 
n(N) /riH ~ 10~ 5 in the pre-shock gas. We also find that P-chemistry is sensitive 
to O- and N-chemistry, because N atoms are destroyed mainly by OH and NO. 
We identify the reactions of O-bearing and N-bearing species that significantly 
affect P chemistry. 

Subject headings: ISMdndividual objects(L1157) — ISM:jets and outflows — shock 
waves — ISM:molecules 

1. Introduction 

Although the solar abundance of phosphorus is 3 x 10~ 7 with respect to H (Asplund et 
al. 2009), only a small number phosphor-containing molecules have been observed in space. 
PN is the only P-bearing species detected in Turner & Bally (1987) detected PN in Ori(KL), 
W51M, and Sgr B2, whereas Matthews, Feldman & Bernath (1987) did not detect PO in 
Ori(KL) or Sgr B2. Considering the non-detection of PO, Turner & Bally (1987) suggested 
that PN is formed by grain disruption. On the other hand, Millar et al. (1987) calculated the 
P-chemistry considering related reactions (e.g., PO + N — > PN + O), and showed that the 
observed PN abundance and non-detection (i.e., the upper limit) of PO can be reproduced 
by gas-phase chemistry at low temperature (50 K). Later, Charnley & Millar (1994) showed 
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that PN becomes abundant within 10 4 yr in hot core models; assuming that P is initially in 
the form of PH 3 ice, they calculated the gas-phase chemistry after the evaporation of PH 3 . 

Recently, Yamaguchi et al. (2011) detected PN in a low-mass star-forming region, L1157 
Bl, for the first time. L1157 is a dark cloud harboring a Class low-mass protostar, which 
drives a well-collimated molecular outflow with the dynamical age of 1.8 x 10 4 yr (Umemoto 
et al. 1992). L1157 Bl is a shocked region formed by interactions between the blue-shifted 
molecular outflow and ambient gases. Since Bl is distant from the protostar, the "pure" 
shock chemistry can be investigated. Many observational studies of the Bl region have been 
conducted to investigate physical and chemical conditions (e.g., Bachiller & Perez Gutierrez 
1997; Hirano & Taniguchi 2001; Sugimura et al. 2011). The kinetic temperature is estimated 
to be 50 — 170 K (Hirano & Taniguchi 2001), and H 2 density is in the range of 10 4 cm -3 
(Hirano & Taniguchi 2001) to 10 5 cm -3 (Bachiller & Perez Gutierrez 1997) in the shocked 
region. Previous work concluded that C-shock occurs at Bl (Gusdorf et al. 2008; Viti et al. 
2011). PN is detected also towards the shocked region L1157 B2. The abundance relative 
to H 2 is estimated to be n(PN)/n(H 2 ) ~ (2 - 6) x 10~ 10 towards Bl and (3 - 7) x 10~ 10 
towards B2 (Yamaguchi et al. 2011). The PN line is blue-shifted, and the line width is as 
broad as 3.8 km s -1 . In addition, PN emission is not detected toward the protostar position. 
Yamaguchi et al. (2011) thus concluded that PN is formed by shock. A subsequent work 
(Yamaguchi et al., in prep) reported that PO is not detected at Bl; the upper limit of 
PO abundance relative to H 2 is 2.5 x 10~ 10 . Detailed modeling of phosphorous chemistry in 
shocks and comparison of the results with these observations will improve our understanding 
of P chemistry. 

In this work, we study the evolution of the P-bearing species in interstellar shocks. 
We solve the chemical reaction network along one-dimensional (ID) C-shock models, which 
mimic the LI 157 Bl region. In §2, we describe the physical and chemical model. Results 
of the numerical calculation are presented in §3. In §4, we discuss the difference between 
our model and that of Charnley & Millar (1994). The dependences of P-chemistry on shock 
parameters are also investigated. We summarize our conclusion in §5. 

2. Model 

We adopt the steady-state ID C-shock model by Jimenez-Serra et al. (2008), in which 
temporal variations of density and temperature (neutral and ion) are given for several pre- 
shock densities and velocities. We make a fitting function of electron temperature, which 
is not included in the model of Jimenez-Serra et al. (2008), based on Draine, Roberge 
& Dalgarno (1983) and Flower et al. (1996). Referring to the observation of L1157 Bl 
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(Bachiller & Perez Gutierrez 1997), we assume the pre-shock density n(H.2) = 1.0 x 10 4 cm -3 
and velocity v — 20 km s _1 . Figure [1] shows the temporal variation of neutral, ion, and 
electron temperatures, H 2 density, and neutral and ion velocities. We define t = as the 
time when the density and temperature start to rise due to the shock. The gas is heated to 
T ~ 1000 K in 10 3 yr, but then cooled to 10 K in 10 4 yr. The gas temperature estimated 
from the observation, 50-170K (Hirano & Taniguchi 2001), is much lower than the peak 
temperature of the shock model because of the beam dilution; the beam size of Hirano & 
Taniguchi (2001) is ~ 0.03 pc, whereas the width of the layer heated to T > 100K is ~ 0.018 
pc in the model. 

Along the flow, we calculate the chemical reaction network, which consists of 658 species 
and 11285 reactions. Chemical reactions and rate coefficients are adopted mainly from the 
OSU network flhttp:// www.physics.ohio-state.edu/~erict); we use the network of Garrod & 
Herbst (2006) at T < 100K and Harada et al. (2010) at T > 100K. We also include collisional 
dissociations in Table Al of Willacy et al. (1998), as well as the reactions of PH3 and PH4 
in Table 1 of Charnley & Millar (1994). Grain surface reactions are included in our chemical 
network to set up the initial abundances of pre-shock matter (see Model B described below). 
In the shock chemistry, however, we neglect the grain-surface reactions, because they are not 
effective in the shocked region within the dynamical time scale of L1157 (1.8xl0 4 yr). 

In C-shock, the neutral and ion temperatures are different. To evaluate the reaction 
rate coefficient in such flows, Flower et al. (1985) introduced effective temperature, T e g, 
including the streaming effect: 

T^boiTeff = ^k ho iT T + ^m in (ui - u n ) 2 

_ mjT n + m n Tj 
rii + m n 

m n mi 

in-m = ■ , 

m n + m,\ 

where u n and u\ are neutral and ion velocities, and m n and m; are neutral and ion masses, 
respectively. We also use the following formula for rate coefficients of collisional dissociation 
of molecules by ions 

k = a x \u\ — u a \ 




if the kinetic energy of the relative motion of the reactants, Et (= \m- m \u\ — u n \ 2 ), is larger 
than the threshold energy E . 



-4- 



In this work, we consider two models. In model A, we assume that heavy-element 
species are completely accreted onto grain surfaces in the pre-shock gas. In model B, we 
calculate molecular evolution in a molecular cloud condition to set the initial abundance of 
the pre-shock gas. 

In hot core models (e.g., Charnley & Millar 1994), it is assumed that heavy elements 
are completely depleted onto grain surfaces initially due to low temperature (10K) and high 
density (n(H 2 ) ~ 10 6 — 10 7 cm -3 ) in prestellar cores. Chemical calculation starts when the 
protostar heats the core and the grain surface species are sublimated. Referring to the hot 
core models, we assumed that heavy elements (heavier than He) are completely depleted on 
grain surfaces in the pre-shock gas in Model A. We adopt the initial molecular abundances 
(Tabled]) of the hot core model of Nomura & Millar (2004). They determined these initial 
abundances to be consistent with the observations of interstellar ices in molecular clouds and 
gas-phase species in hot cores. Previous studies (e.g., Draine, Roberge & Dalgarno 1983; 
Caselli, Hartquist & Havnes 1997) showed that the sputtering of grain-surface species is 
effective in C-shock with 20 km s -1 . Jimenez-Serra et al. (2008) calculated the sputtering 
of dust surface species and showed that the time scale of sputtering is about 100 yr, if the 
initial H 2 density is 10 4 cm -3 . Since this sputtering time scale is much shorter than the 
heating time scale of neutral gas and electrons (Figure [1]) and the chemical time scale (see 
§3), we simply start our calculation with all the molecules desorbed to the gas phase at t =0. 

In model B, we set the initial molecular abundances by calculating the chemical evo- 
lution in a molecular cloud condition (n(H 2 )=10 4 cm -3 , T = 10 K), which is the same as 
the parameters in our pre-shock gas. Elemental abundance is the same as Model A. The 
molecular cloud calculation starts with atoms and ions except for hydrogen and phosphorus; 
H is in molecular form (Table |2]) and P is in the solid phase (see below). We adopt the 
molecular abundances at t = 1 x 10 5 yr as the initial condition of the shock chemistry. At 
the start of the shock chemistry, all ice-mantle species formed in the molecular cloud are 
desorbed to the gas phase. 

The initial form and abundance of phosphorous relative to hydrogen require some ex- 
planation. Although the solar abundance of phosphorous is 3 x 10~ 7 , P-bearing species are 
not detected in quiescent molecular clouds (e.g., Turner et al. 1990). If all P is in the 
gas phase (i.e., P + ) initially, theoretical models predict that PN abundance becomes high 
enough to be easily detected in molecular clouds. A significant amount of phosphorous thus 
should be depleted in the solid phase in quiescent clouds (e.g., Millar et al. 1987), but the 
dominant form of P in solid is unknown. Charnley & Millar (1994) assumed that PH 3 ice is 
as abundant as 1.2 x 10~ 8 in the prestellar core, and that only this amount of phosphorous 
is desorbed to the gas phase upon heating via star formation. Following Charnley & Millar 
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(1994), we assume PH 3 abundance of 1.2 x 10 8 in Model A. 

In Model B, we assume that all P is in the solid phase in molecular clouds. We then 
assume that a fraction of P is desorbed to the gas phase in atomic form via sputtering in the 
shock. To make the comparison with Model A easy, the desorbed abundance of P atoms is 
set to 1.2 xlO -8 . We note that the abundances of P-bearing species are proportional to the 
assumed elemental abundance of P (P atom or PH 3 initially) in the shocked gas. Although 
we assume a gaseous P abundance of 1.2 x 10~ 8 in our fiducial model, we vary this abundance 
in order to reproduce the observed abundance of PN and the non detection of PO. The result 
tells us how large of a fraction of P is in the relatively volatile component in the solids that 
is desorbed to the gas phase via sputtering in the C shock. 



3. Result 
3.1. Model A 

Figure [2] shows the time evolution of P-bearing species in Model A. Abundances are 
relative to the number density of hydrogen nuclei (%). Major chemical reactions of P- 
bearing species are shown in Figure [3j At first, P is in PH3, which reacts with H atoms, 
converting them to P atoms: 

H + PH 3 ^PH 2 + H 2 (1) 

H + PH 2 -> PH + H 2 (2) 

H + PH ->• P + H 2 . (3) 

These reactions have activation barriers of 735K, 318K, and 416K, respectively, and thus 
become active only when the gas temperature is higher than ~ 100 K, i.e., from t =several 
10 2 yr to a few 10 3 yr in our model. P atoms react slowly with H 3 + to form PO 

P + H 3 + -> HPO+ + H 2 (4) 

HPO + + H 2 O^PO + H 3 + . (5) 
PN is then formed by the following reactions 

PO + N-^PN + O (6) 

PH + N^PN + H. (7) 



For comparison with the observation, we refer to the abundances at the dynamical time 
scale of the outflow of L1157, ~ 10 4 yr (Umemoto et al. 1992). In Figure [21 PN abundance 
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is lower, and PO abundance is much higher than observed in LI 157 Bl. Considering the 
reactions ([6]) and ([TJ) , it is clear that the abundances of PN and PO depend on that of N 
atoms. In the current model, the dominant form of nitrogen is N 2 , and N atoms are formed 
by the collisional dissociation of N 2 

N 2 + Fe+ ->■ N + N + Fe + 

and the following reactions 

He + + N 2 -»■ N + + N + He 

N 2 + cosmic ray — > N + N. 

Collisional dissociation is effective only at t ~ 200 yr. The latter two reactions proceed 
slowly and thus produce insufficient amounts of N atoms to form PN. We note here that the 
initial zero abundance of N atoms would be artificial. According to theoretical models and 
observations (Aikawa et al. 2001; Maret et al. 2006), nitrogen stays in an atomic form for 
a significant time scale (> 10 5 yr) in molecular clouds. Therefore, we next investigate the 
dependence of PN and PO abundances on the initial abundance of N atoms. 

Figure H] shows the temporal variation of PN, PO, and N atoms with various initial N 
atom abundances. The elemental abundance of nitrogen is kept constant; in a model with 
high N atom abundance, the initial N 2 abundance is low. As expected, PN abundance is 
significantly increased, whereas PO abundance is lowered, if the initial N atom abundance 
is high. The abundances of other P-bearing species (P atom, PH, and PH 3 ) do not change 
significantly from those shown in Figure [2j N atom decreases at t = 10 4 — 10 5 yr (Figure H] 
c); it is converted to N 2 by the following reactions 

N + OH NO + H (8) 

N + NO -> N 2 + O (9) 

NO + NH 2 -»■ N 2 + H 2 0. 

Then PN is converted back to PO as N atom decreases (Figure [3]). 

If the initial N atom abundance is 1.0 x 10~ 6 , the observed abundance of PN is reproduced 
on a relatively limited time scale ~ 10 4 yr, but PO abundance exceeds the observed upper 
limit. On the other hand, if the initial N atom abundance is 1.0 x 10~ 5 , PO abundance is 
consistent with the upper limit for a longer time scale, but PN abundance slightly exceeds 
the observed abundance ((1 — 3) x 10 -10 relative to n H ) at t ~ 10 4 yr. We note here that 
these abundances of P-bearing species are proportional to the initial PH3 abundance. In the 
current calculation, we assumed the initial PH3 abundance of 1.2 x 10~ 8 following Charnley 
& Millar (1994). However, this assumption is rather arbitrary; Charnley & Millar (1994) did 
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not compare their PN abundance with observation quantitatively, because the abundance 
estimates in hot cores were uncertain due to beam dilution and difficulty in estimating the H 2 
column density. Thus far, there are no direct observational constraints on PH 3 abundance. 
Dot-dashed lines in Figure H] show the temporal variation of PO and PN with the initial PH 3 
abundance of 1.2 x 10~ 9 and N atom abundance of 10~ 5 ; both PN and PO abundances are 
consistent with the observation. 

We also note that the initial form of desorbed phosphorus is not necessarily PH 3 , since 
there are no observational constraints on P in solids. Our results do not change if P is 
desorbed as P atoms, instead of PH 3 , because P atoms are on the path to PN in the chemical 
network (see Figure |3]). 

3.2. Model B 

In Model A and Figure El we can see that N atoms, O atoms, and protonated water 
are major reactants with P-bearing species. This indicates that P chemistry may depend 
on the initial abundances of water and O atoms, as well as on the N atom abundance. In 
order to check the dependence of P chemistry on initial molecular abundances, here we set 
the initial abundances by calculating molecular evolution in a molecular cloud condition, 
n(H 2 ) = 1.0 x 10 4 cm -3 and T = 10 K (Figure E]). We adopt the molecular abundance at 
1 x 10 5 yr as the initial abundance for the shock model. Ice mantle species are assumed to 
be desorbed to the gas phase as in Model A. We assume that phosphorus is totally depleted 
in the solid phase in molecular clouds, but a fraction (1.2 x 10~ 8 relative to hydrogen nuclei) 
of P is desorbed to the gas phase in atomic form via sputtering in the shock. 

Figure [6] (a) shows the temporal variation of P-bearing species in the C shock model. 
Major chemical reactions of P-bearing species are mostly the same as in Model A. Since 
N atoms are as abundant as 10 -5 initially, PN becomes abundant, whereas PO is depleted 
at t ~ 10 4 yr. If we adopt the abundances at t — 1 x 10 6 yr in the molecular cloud as 
initial conditions for the shock chemistry, PO abundance exceeds the observed upper limit 
at t ~ 10 4 yr (Figure [6b), because the initial N atom abundance is low (~ 7.0 x 10~ 7 ). 

As in Model A, the abundances of P-bearing species are proportional to the assumed 
elemental abundance of P in the gas phase. If we adopt the molecular cloud abundance at 
10 5 yr and the initial P atom abundance of ~ 3 x 10~ 9 , our PN and PO abundances become 
consistent with the observation at ~ 10 4 yr. 
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4. Discussion 

4.1. Differences between our results and those of Charnley & Millar (1994) 

In our model, PN becomes abundant only if initial N atom abundance is high (> 10~ 5 ). 
In the hot core model of Charnley Sz Millar (1994), on the other hand, abundant PN is formed 
in ~ 10 4 yr even though N is all in N2 initially. Why do we need abundant N atoms? Is it 
due to the difference in physical models (hot core vs. shock) or chemical reaction networks? 
We note that Charnley & Millar (1994) used the network model of UMIST91 (Millar et al. 
1991). Since reaction rate coefficients have been updated significantly since UMIST91, it is 
worth investigating how our network model differs from UMIST91, and if our conclusion, 
that there is high N atom abundance in pre-shock gas, sensitively depends on the reaction 
network model. 



4.I.L Comparison between chemical network models 

In order to compare our network model with UMIST91, we calculated the same hot 
core model as Charnley & Millar (1994). We adopted the initial abundance in their Table 2; 
it is essentially similar to our initial condition of Model A with abundant N 2 and saturated 
molecules such as H 2 0. Temperature and density are set to T=100 K and n(H 2 ) = 10 7 cm -3 . 

The result is shown in Figure [7] (a). At t ~ 10 4 — 10 5 yr, PN abundance is much lower, 
and PO abundance is much higher than that in Charnley & Millar (1994), in which PO 
decreases and PN abundance reaches several 10 -9 at t ~ 10 4 yr (see Figure 1 of Charnley 
& Millar 1994). Nitrogen atoms, which are needed to form PN, are not abundant until t ~ 
several 10 5 yr; although N atoms are formed by N 2 + He + , their abundance is suppressed by 
reactions (|Sj) and Q. Although N atom abundance is not shown in Charnely & Millar (1994), 
it might be higher than in our model. Since N atom abundance depends on N-chemistry and 
O-chemistry (see reactions (jHD and we compared the rate coefficients of major N- and 
O- reactions. Table 3 lists several reactions that we found to have different rate coefficients 
in our network and UMIST91. For example, the rate coefficient of 

CO + OH^C0 2 + H (10) 

is 500 times higher in UMIST91 than in our model at T=100K; the rate coefficient in 
UMIST91 is adopted from Smith (1998), whereas our rate coefficient is from Frost et al. 
(1993). For reactions that are not included in UMIST91, we put in the column of A 
parameter of "UMIST91" in Table [3J It is natural that most of the reactions listed in 
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Table [3] are neutral-neutral reactions, whose rate coefficients are hard to measure at low 
temperatures in laboratories. 

We adopted the rate coefficients of UMIST91 for the reactions in Table 3 (hereafter, 
modified network model) and calculated the hot core chemistry of Charnley & Millar (1994) 
(Figure [7b). Now, the temporal variation of P-bearing species is more similar to that in 
Charnley & Millar (1994) than to that in Figure [TJ (a). We can conclude that there are 
differences in N-chemistry and O-chemistry between our network and UMIST91, which sig- 
nificantly affect P-chemistry. 

4-1.2. Shock chemistry with the modified network model 

Using the modified network model, we calculated the C-shock model to see how the 
evolution of P-bearing species differs from our original model. The initial abundances are 
the same as in Model A. Figure E] shows the PN and PO abundances with various initial 
N atom abundances (0,1.0 x 10~ 7 , 1.0 x 10~ 6 , and 1.0 x 10~ 5 ). Although these results 
are slightly different from those in Figure HI our original conclusion stands: PN and PO 
abundances become consistent with the observed value and the upper limit only if N atom 
abundance is as abundant as > 10~ 5 . 



4-1.3. Dependence of hot core model on physical parameters 

Using the modified network model, PN becomes abundant at ~ 10 4 yr in the hot core 
model (as in Charnley & Millar 1994) even if the initial N atom abundance is zero. In the 
shock model, on the other hand, we still need initially high N atom abundance to make PN 
abundant at ~ 10 4 yr. Since we are now comparing the two models with the same chemical 
reaction network, the difference should be attributable to the physical conditions. In the 
C-shock model for L1157 Bl, density (^(H^) = 10 4 cm -3 ) is much lower than that in the 
hot core model (^(H^) = 10 7 cm -3 ), and the temperature varies temporally from 10 K to 
1000 K, whereas the hot core model has a fixed temperature (100 K). In the following we 
calculate hot core models with various densities and temperatures in order to determine 
density, temperature, or both affect the P chemistry. We use the modified reaction network 
and the initial conditions shown in Table 2 in Charnley & Millar (1994). 

Figure E(a) shows the PN abundance in hot core models of T — 100 K with the gas 
densities of n(H 2 ) = 1.0 x 10 4 cm -3 , 1.0 x 10 5 cm -3 , and 1.0 x 10 7 cm -3 . We can see 
that PN is formed in shorter time scales at higher densities. In general, chemical reactions 
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proceed faster at higher densities. Whereas the collisional time scale of gaseous particles 
is proportional to n -1 , the chemical time scales dependence on density is weaker than n^ 1 , 
because the ion-molecule reactions play major roles in the reaction network and because the 
ionization degree is proportional to n" 1 / 2 . 

Figured] (b) shows the PN abundance in hot core models of n(B.2) = 1-0 x 10 7 cm -3 with 
T = 10 K, 30 K, 50 K, and 100 K. At lower temperatures, PN is less abundant, because the 
N atom is destroyed by the reaction with the OH (reaction®). Figure M^c) shows the OH 
abundance in the same models as in Figure[H](b). OH increases with time at low temperatures 
(10 K and 30 K), whereas at T — 100 K, OH is destroyed by the reaction (ITU1) . which has 
an activation barrier of 390 K in UMIST91. We note that in the C-shock model (Figure 1), 
temperature initially increases to ~ 1000 K, but decreases to ~ 10K after several 10 3 yrs, 
when the OH abundance increases and hampers the formation of PN. 

From these discussions, we can conclude that both high density and high (constant) 
temperature help to make PN abundant in the hot core model, whereas in the C-shock 
model, we need high N atom abundance in the pre-shock gas in order to make PN abundant 
at 10 4 yr. 



4.2. Dependence on shock parameters 

So far, we have adopted the C-shock model of n(H2) = 1.0xl0 4 cm -3 and v = 20 km 
s _1 (Figured]). In this section, we investigate the dependence of P-chemistry on C-shock 
parameters: pre-shock density and velocity. We use our original reaction network and adopt 
the initial conditions in Model A. The initial N atom abundance is 1.0 xlO -5 . Figure [TD] 
shows the temporal variation of neutral, ion, and electron temperatures, H2 density, and 
neutral and ion velocities in two C-shock models: (a, c) the pre-shock density n(H2) = 10 4 
cm -3 and velocity v — 40 km s _1 , and (b, d) n(H 2 ) = 10 5 cm -3 and v — 20 km s^ 1 . The 
peak temperature is higher in the model with higher velocity, whereas the cooling time is 
shorter in the model with higher pre-shock density. 

Figure [TT1 (a) shows the temporal variation of P-bearing species in the model of n(H.2) = 
10 4 cm -3 and v — 40 km s _1 . Since the peak temperature is higher than the model with 
v = 20 km s _1 , N atoms are depleted in the high-temperature region (T ~ 2000K) by 

H 2 + N — )• NH + H, (11) 

which has an energy barrier of 1.66 x 10 4 K. Then, PN is less abundant and PO is more 



- 11 - 



abundant than Figure |4] (b). NH is further hydrogenated to NH3 by 



NH + H 2 ->■ NH 2 + H 



(12) 



NH 2 + H 2 -> NH 3 + H, 



(13) 



which also have activation barriers. At high temperatures, H atoms increase due to collisional 
dissociation of H 2 , which starts the back reactions of fill I) - 0131) . However, back reactions are 
not efficient enough to stop the decline of N atom abundance. Figure [TT1 (b) shows the same 
shock model but with the initial N atom abundance of 7.46 xlO -5 ; nitrogen is initially all 
in atomic form. PN is more abundant than PO at ~ 10 4 yr. If we reduce the initial PH3 
abundance by a factor of a few, PO abundance can be marginally lower than the upper limit, 
whereas PN abundance is in the range of observed values. 

Figure [11] (c) shows, on the other hand, the temporal variation of P-bearing species 
with the pre-shock density n(H 2 ) = 1.0xl0 5 cm~ 3 and velocity v = 20 km s _1 . The results 
are similar to those in Figure HJ PN is abundant and PO is depleted at t ~ 10 4 yr. 

In summary, the observed abundances of PN and PO are reproduced in the C-shock 
models if the peak temperature is less than T < 2000 K. In the shock model with higher 
peak temperatures, only marginal agreement is found with an extremely high initial N atom 
abundance. 



In this work, we have calculated the phosphorus chemistry in ID steady-state C-shock 
models to examine how PN is formed and why PO is not detected in the shocked region 
LI 157 Bl. We found that PN and PO abundances in our model become consistent with the 
observation, if the following conditions are satisfied: 

• a fraction of phosphorus, ~ 10~ 9 relative to hydrogen, is desorbed from the grain to 
the gas phase as PH 3 or P atom via sputtering in the shock 

• initial N atom abundance is high (> 10~ 5 ) 

• the peak temperature of the shocked gas is < 2000 K. 



5. 



Conclusions 



We have also compared our reaction network model with UMIST91, which is adopted in 
the hot core model of Charnley & Millar (1994). Phosphorus chemistry sensitively depends 
on the abundance of N atoms, which in turn depends on N- and O-chemistry. We identified 
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several neutral-neutral reactions with different rate coefficients in our network models and 
UMIST91, and these reactions significantly affect the P-chemistry. We confirmed that our 
conclusion — a high abundance of N atoms in the pre-shock gas is needed to make PN 
abundant — is valid even if we modify the rate coefficients of these reactions to those of 
UMIST91. Furthermore, our results indicate that, even in hot core models, high N atom 
abundance is needed to make PO less abundant than PN at t ~ 10 4 yr, if we adopt our 
updated network instead of UMIST9f . 

We are grateful to NRO survey team, especially to S.Yamamoto and T.Yamaguchi 
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Species 


Abundance 


Species 


Abundance 


H 2 


0.5 


H 2 S 


l.Ox 10~ 7 


He 


0.1 


OCS 


5.0x 10~ 8 


H 2 


2.8x 10- 4 


Si 


3.6x 10- 8 


CO 


1.3x 10- 4 


PH 3 


1.2x 10~ 8 


H 


5.0x 10~ 5 


S 


5.0x 10~ 9 


N 2 


3.7x 10~ 5 


C 2 H4 


5.0x 10~ 9 


co 2 


3. Ox 10~ 6 


C 2 H 5 OH 


5. Ox 10~ 9 


H 2 CO 


2.0x 10~ 6 


C 2 H 6 


5.0x 10~ 9 


o 2 


l.Ox 10~ 6 


Fe+ 


2.4x 10~ 8 


NH 3 


6.0x 10~ 7 


m 


l.Ox 10" 9 


C 2 H 2 


5.0x 10~ 7 


H+ 


l.Ox lO" 11 


CH 4 


2.0x 10- 7 


He+ 


2.5x 10~ 12 


CH3OH 


2.0x 10- 7 







Table 1: Initial molecular abundances relative to hydrogen nuclei in Model A 



Species Abundance 


Species Abundance 


H 2 5.0x 10- 1 
He l.Ox 10- 1 

4.18x 10~ 4 
C+ 1.36x 10" 4 

N 7.46 x 10" 5 


H 5.0x 10" 5 
S+ 1.55x 10- 7 
Si+ 3.6xl0- 8 
Fe+ 2.4x 10~ 8 



Table 2: Initial abundances in the molecular cloud model for Model B 
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UMIST91 


reaction 
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k(T — 1 (10 K1 a 

l\i \ -L lUu 1 V J 


A 




Q 


k(T — inn to 


CO + OH -> CO2 + H 


2.81 (-13) b 





176 


4.8 (-14) 


3.1 (-10) 


-1.15 


390 


2.2 (-11) 


O 4- NHo -4 NO 4- Ho 


l.U (-11J 


u 


u 


l.U (,-llJ 


u 








O 4- NH-i — ► NH 4- OH 


2.0 (-11) 








2.0 (-11) 


3.5 (-12) 


0.5 





9 09 ( 1 9^ 


O + NH -4 NO 4- H 


1.16 (-10) 








1.16 (-10) 


1.73 (-11) 


0.5 





1 Or 1 "M 
l.U ^-iij 


O 4- PH -4 PO + H 


1.0 (-10) 








1.0 (-10) 


4 C-1 1 1 

l.U ^ 11) 








4 D f- 1 1 1 


HPO+ 4- H2O — > H 3 0+ + PO 


_c 


_c 


_c 




1 i Q\ 
l.U 








1 ( Ql 
i.u ^-yj 


N + NO — > N 2 + O 


3.0 (-11) 


-0.6 





5.8 (-11) 










^ ZL T 11 s ! 

0.^ v -ii; 


N + OH -4 NO + H 


7 5 C-1 1 1 


-0.18 





9.14 (-11) 


5.0 (-11) 


0.5 





9 SQ C-1 1 ^ 

i.Oil 1 111 


N 4- CH 3 -> H 2 CN + H 


8.6 (-11) 








8.6 (-11) 











O + OH — > O2 4- H 


3.5 (-11) 








3.5 (-11) 


7 9 (-11) 








7 9 (-11 s ) 


O + HNO -> OH + NO 


3.8 (-11) 








3.8 (-11) 


1.44 (-11) 


0.5 





8.31 (-12) 


O 4- NH 2 -4 HNO + H 


8.0 (-11) 








8.0 (-11) 


3.2 (-12) 








3.2 (-12) 


CH5O+ + c -4 CH 3 + OH + H 


4.64 (-7) 


-0.67 





9.69 (-7) 











CH s O+ + c -4 CH 3 + H 2 O 


8.19 (-8) 


-0.67 





1.71 (-7) 











CH 4- H 2 -4 CH 3 


3.25 (-17) 


-0.6 





6.28 (-17) 











CH3OCH+ + c -4 CH 3 4- CH 4 4- O 


1.5 (-7) 


-0.5 





2.6 (-7) 











CH3OCH+ + e -4 CH4O 4- CH 3 


1.5 (-7) 


-0.5 





2.6 (-7) 











O + C3H3 -4 C 2 H 2 + HCO 


2.31 (-10) 








2.31 (-10) 











O + CH 2 -4 HCO + H 


5.01 (-11) 








5.01 (-11) 


1.44 (-11) 


0.5 


2000 


1.71 (-20) 


H + HCO -4 H 2 + CO 


1.5 (-10) 








1.5 (-10) 


3.0 (-10) 








3.0 (-10) 


PO + N -4 PN + O 


3.0 (-11) 


-0.6 





5.8 (-11) 


3.4 (-11) 








3.4 (-11) 



Table 3: List of reactions with different rate coefficient in our model and UMIST91 which 
affect the P-chemistry. Rate coefficients are given as a function of A,B, and C: k = 
A(T/3W) B exp(-C/T). 



"Rate coefficient at T = 100 K. 
b a(-b) means a x 10~ b . 

c The rate coefficient in OSU model is calculated using equations (3), (5) - 



(7) of Harada et al. (2010) 
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Fig. 1. — (a) Temporal variation of neutral, ion, and electron temperatures and H2 density 
in the C-shock model with pre-shock density n(H2) = 1.0 x 10 4 cm -3 and velocity v = 20 km 
s -1 . (b) Temporal variation of the neutral and ion velocities in the C-shock model. Velocity 
is in the frame co- moving with the shock front. Dot-dashed line depicts the ion-neutral drift 
speed (\m —u n \). 
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Fig. 2. — Temporal variation of P-bearing species in Model A. Initial N atom abundance is 
set to zero. 
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Fig. 3. — Reaction network of P-bearing species in the shocked gas. 
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Fig. 4. — PN, PO, and N atom abundances in Model A with various initial N atom abun- 
dances, which are each attached to a line. Dot-dashed lines represent PN and PO abundances 
with the initial Phh abundance of 1.2xl0 -9 and N atom abundance of 10 -5 . 
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Fig. 5. — Molecular evolution in the molecular cloud condition, nn = 2.0 x 10 4 cm 3 and 
T = 10 K. 
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n(i)/n H ( a ) n(i)/n H (b) 




Fig. 6. — Temporal variation of P-bearing species in Model B. We adopt the molecular cloud 
abundances at (a) t = 1 X 10 5 yr and (b) t = 1 X 10 6 yr as the initial conditions of shock 
chemistry. P atoms are assumed to be desorbed to the gas phase by sputtering in the shock. 



n(i)/n H n(i)/n H 




time(yr) time(yr) 

Fig. 7. — Temporal variation of P-bearing species in the hot core model of T=100K and 
nn=2.0xl0 7 cm -3 . We use our original network model in (a) and the modified network 
model in (b). 
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Fig. 8. — PN and PO abundances in the shock model as shown in figure H] but with the 
modified network model. 
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Fig. 9. — (a) Dependence of PN abundance on gas density (n H =2.0 x 10 7 , 2.0 x 10 5 , 2.0 x 10 4 
cm -3 ), and (b, c) dependence of PN and OH abundances on temperatures (T=10, 30, 50, 
100K) in hot core model. 
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Fig. 10. — Temporal variation of temperature and density in the C-shock model with (a) pre- 
shock density ^(H^) = 1.0 x 10 4 cm -3 and velocity v = 40 km s _1 , and (b) ra(H2) = 1.0 x 10 5 
cm -3 and velocity v — 20 km s _1 . (c-d) Temporal variation of the neutral and ion velocities 
in the two C-shock models. 
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Fig. 11. — Temporal variation of P-bearing species in the C-shock model with pre-shock 
density nn = 2.0 x 10 4 cm -3 and velocity v = 40 km s _1 (a-b), and % = 2.0 x 10 5 cm -3 
v = 20 km s _1 (c). Initial N atom abundance is 1.0 x 10~ 5 in (a) and (c), whereas it 
is 7.46 x 10 -5 in (b). Dashed and dot-dashed lines in panel (b) depict the observed PN 
abundance and upper limits on PO abundance, respectively. 



